Let us join dialogues with information contained in characters and episodes. This allows us to know the gender of the most talkative characters and the season of each episode. Notice that dialogues terminates at episode 568 (episode 16 of season 26), whereas episodes at episode 596 (last episode of season 27).

dialogues <- dialogues %>%
  left_join(characters, by = c("role" = "name")) %>%
  mutate(gender = fct_explicit_na(gender, na_level = "Unknown")) %>%
  left_join(episodes[,c("season", "episode_id")], by = "episode_id")

dialogues <- dialogues %>%
  mutate(season = ifelse(is.na(season), 8, season)) %>%
  select(line_number, season, number = episode_id, location, role, gender, line)

Let us also join episodes with information from guests. This will be useful to get, for instance, the guest star names for each episode (if any), and the roles they played. Note that guests terminates at episode 662 (episode 23 of season 30).

episodes <- episodes %>%
  left_join(guests, by = c("number", "season", "prod_code", "title")) %>%
  mutate(guest = ifelse(!is.na(role), TRUE, FALSE))

Lastly, we tidy dialogues, which conveniently puts the speaking lines into a one-word-per-row format.

dialogues.tidy <- dialogues %>%
  unnest_tokens(word, line)

The Simpsons is characterized by a marked gender imbalance, which is also reflected in the show’s writing staff. More than 75% of the characters with recorded gender are male. The only female leading characters are Marge and Lisa Simpson. In contrast, among the supporting cast, we find Edna Krappabel-Flanders (the teacher at Springfield Elementary School), and the twins Selma and Patty Bouvier (Marge’s older sisters).

lollipop <- function(df, x, fill.var, count = TRUE, index = NULL, title = "", ylab = "",
                     sub = NULL, pos = "bottom", labels = TRUE){

  x        <- enquo(x)
  fill.var <- enquo(fill.var)

  count_function <- function(df, x, y){

  x <- enquo(x)
  y <- enquo(y)

  df %>%
    count(!! x, !! y, sort = TRUE)
  }

  if(count){
    df <- df %>%
      count_function(x = !! x, y = !! fill.var)
  }

  if(!is.null(index)){
    df <- df %>%
      slice(index)
  }

  p <- df %>%
    ggplot(aes(!! x, n, fill = !! fill.var)) +
    geom_segment(aes(
      x = fct_reorder(!! x, n),
      xend = !! x,
      y = 0,
      yend = n,
      color = !! fill.var),
      size = 0.8,
      show.legend = FALSE) +
    geom_point(aes(
      color = !! fill.var),
      size = 4,
      alpha = 0.6) +
    coord_flip() +
    labs(x = "", y = ylab, title = title, subtitle = sub) +
    scale_fill_manual(name  = firstup(as_label(fill.var)),  values = rev(colors[1:2])) +
    scale_color_manual(name = firstup(as_label(fill.var)),  values = rev(colors[1:2])) +
    theme(legend.position = pos)

  if(labels){
    p <- p +
      scale_y_continuous(labels = comma)
  }
  p
}

plot.top.chars <- dialogues.tidy %>%
  lollipop(x  = role, fill.var = gender, index = c(1:10), title = "Top 10 Simpson leading characters",
           ylab = "Number of spoken words in 26 seasons")
plot.supporting.chars <- dialogues.tidy %>%
  lollipop(x  = role, fill.var = gender, index = c(11:35), title = "The Simpsons supporting cast",
           ylab = "Number of spoken words in 26 seasons")

ggarrange(plot.top.chars, plot.supporting.chars, nrow = 1, common.legend = TRUE, legend="bottom",
          widths = c(0.95, 1.05))

Let us look at how the total number of words of the top four leading characters (Homer, Marge, Bart, and Lisa Simpson) has evolved throughout the seasons.

Top20Chars <- dialogues.tidy %>%
  count(role, sort = TRUE) %>%
  distinct(role) %>%
  head(20) %>%
  pull(role)

dialogues.tidy %>%
  filter(role %in% Top20Chars[1:4]) %>%
  group_by(season, role) %>%
  summarize(nwords = n()) %>%
  ungroup() %>%
  mutate(role = reorder(role, -nwords)) %>%
  ggplot(aes(x = season, y = nwords, fill = role, color = role)) +
  scale_y_continuous(labels = comma) +
  scale_x_continuous(breaks = seq(1, 26, 5)) +
  geom_col(alpha = 0.65) +
  facet_wrap(~ role, nrow = 1) +
  labs(x = "Season", y = "Number of spoken words per season",
       title = "Homer Simpson has dominated all dialogues in every season") +
  scale_fill_manual(name = "Character", values = palette[c(8,14,6,1)]) +
  scale_color_manual(name = "Character", values = palette[c(8,14,6,1)]) +
  theme(legend.position = "bottom")

We notice an overall downward trend for all three indicators. Older seasons seem to be the most appreciated, most rated, and most viewed. As a matter of fact, The Simpsons show received acclaim throughout its first nine or ten seasons, which are generally considered its “Golden Age”, but has been criticized for a perceived decline in quality over the years.

To be fair, the availability in recent years of a variety of network channels and Internet streaming platforms has caused a systematic drop in the number of TV views not only for The Simpsons but also for many other shows.

Let us look at the pairwise-relationships between IMDb ratings, IMDb votes, and TV views. There seem to be positive linear relationships between each pair, suggesting that when one of those indicators increases, so do the other two.

line_smooth <- function(data, mapping, method="loess", ...){
  p <- ggplot(data = data, mapping = mapping) +
    geom_point(...) +
    geom_smooth(method=method)
  p
}
bi_density_plot <- function(data, mapping, palette = 4, ...){
  p <- ggplot(data, mapping = mapping) +
    stat_density_2d(aes(fill = stat(level)), geom = "polygon") +
    scale_fill_distiller(palette = palette, direction = 1)
  p
}

episodes %>%
  select("IMDb rating" = rating, "IMDb votes" = votes, "TV views in US (millions)" = us_views) %>%
  ggpairs(upper = list(continuous = wrap(bi_density_plot, palette = 1, size = 0.2)),
          diag  = list(continuous = wrap("barDiag", fill = "#8d99ae", bins = 27)),
          lower = list(continuous = wrap(line_smooth, color = "black", size = 0.2))) +
  ggtitle(label = "Pairwise variable comparisons")

The correlation matrix of episodes certainly shows strong inverse relationships between the episode number and the ratings, the votes, and the views. This indicates that more recent seasons generally have lower values of those performance indicators.

episodes %>%
  select("Episode number" = episode_id, "IMDb rating" = rating, "IMDb votes" = votes, "TV views in US" = us_views) %>%
  cor() %>%
  ggcorrplot(type = "lower", colors = c("#6D9EC1", "white", "#E46726"), outline.col = "white",
             legend.title = "Correlation", lab = TRUE, ggtheme = ggplot2::theme_light()) +
  labs(title = "Correlation matrix")

Who are the guest stars who played themselves in multiple episodes? At the top, we find the physicist and cosmologist Stephen Hawking with 4 appearances across 30 seasons, followed by the comic-book writer Stan Lee, the filmmaker Ken Burns, and the actor Gary Coleman, all of them with 3 occurrences. The gender imbalance in the original characters is also reflected in the guest appearances, with just 5 women playing themselves twice in 30 seasons. The majority of the guest stars, however, just appears in a single episode.

guests %>%
  filter(self) %>%
  count(guest_star, sort = TRUE) %>%
  filter(n > 1) %>%
  add_column(gender = c(rep("Male", 7), "Female", rep("Male", 7), rep("Female", 2), "Male",
                        "Female", rep("Male", 14), rep("Female", 1), rep("Male", 2))) %>%
  lollipop(x = guest_star, fill.var = gender, ylab = "Number of guest appearences in 30 seasons",
           count = FALSE, title = "Who has played themselves in multiple Simpsons episodes?")

It would be interesting to see whether the number of guest appearances has changed over time. In the show’s early years, most guest stars have voiced original characters, but as the show has continued, the number of those appearing as themselves has increased, especially throughout seasons 12 to 19. In more recent seasons, the gap between the two conditions seems to have become more pronounced.

guests %>%
  mutate(self = factor(self, levels = c(FALSE, TRUE), labels=c("Playing an original character", "Playing themselves"))) %>%
  group_by(season, self) %>%
  summarize(n = n()) %>%
  ggplot(aes(season, n, color = self)) +
  geom_line(size = 1.2) +
  scale_color_manual(name = "Guest star", values = colors[3:4]) +
  scale_x_continuous(breaks = seq(1, 30, 2)) +
  theme(legend.position = "bottom") +
  labs(x = "Season", y = "Number of guest appearences",
       title = "The number of guest appearances over seasons")

When guest stars are voicing a character, how long do they talk? Let us explore the average number of lines reserved for guest stars. The guest stars with the most lines per episode are the ones voicing a narrator or an announcer and are usually not playing themselves. The only exception of a guest star playing themselves and having a fair amount of lines has been Lady Gaga.

guests <- guests %>%
  mutate(role = ifelse(self, guest_star, role))

guests_summarized <- guests %>%
  filter(season <= 27) %>%
  group_by(guest_star, role, self) %>%
  summarize(nb_episodes  = n())

guest_roles <- guests_summarized %>%
  inner_join(dialogues %>%
               count(role, sort = TRUE, name = "nb_lines"),
             by = "role") %>%
  mutate(lines_per_episode = nb_lines/ nb_episodes)

guest_roles %>%
  arrange(desc(lines_per_episode)) %>%
  show_table(head = 15, align = c("l", "l", rep("r", 4)),
             col.names = c("Guest star", "Role", "Playing themselves", "Number of episodes",
                           "Number of lines", "Number of lines per episode"),
             caption = "Top 15 guest stars with the highest number of lines per episode")

Guest stars playing themselves tend to have fewer lines than those playing an actual character on the show.

guest_roles %>%
  mutate(self = ifelse(self, "Playing themselves", "Playing an original character")) %>%
  ggplot(aes(lines_per_episode)) +
  geom_histogram(aes(fill = self), binwidth = 2, center = 1, show.legend = FALSE) +
  facet_wrap(~ self, ncol = 2) +
  scale_fill_manual(values = colors[3:4]) +
  labs(x = "Number of lines per episode", y = "Frequency",
       subtitle = "Most guest stars, especially those playing themselves, have relatively few lines per episode")

Is there a difference in IMDb ratings, IMDb votes, and TV views in the episodes with guest stars playing themselves versus playing an original character? Somehow. It seems that episodes with guests starring themselves have, on average, lower ratings, votes, and views than the episodes with no guest star. According to the two-samples Wilcoxon test, the differences in the mean levels are statistically significant at a 5% level.

episodes_by_self <- episodes %>%
  filter(!is.na(self)) %>%
  group_by(self)

episodes_by_self %>%
  summarize(avg_rating = mean(rating),
            avg_votes  = mean(votes),
            avg_views  = mean(us_views)) %>%
  mutate(self = ifelse(self, "Playing themselves", "Playing an original character")) %>%
  add_row(self = "p-value",
          avg_rating = wilcox.test(rating ~ self, data = episodes_by_self, exact = FALSE)$p.value,
          avg_votes = wilcox.test(votes ~ self, data = episodes_by_self, exact = FALSE)$p.value,
          avg_views = wilcox.test(views ~ self, data = episodes_by_self, exact = FALSE)$p.value) %>%
  show_table(col.names = c("Guest star", "IMDb rating", "IMDb votes", "TV views in US (millions)"),
             align = c("l", rep("r", 3)), digits = 3,
             caption = "Average performance indicators and p-values of Wilcoxon rank sum test")

plot_violin <- function(df, x, y, ylab, title = "", limits = NULL,
                        breaks = NULL, labels = comma){

  x <- enquo(x)
  y <- enquo(y)

  data_summary <- function(x) {
    m <- mean(x)
    ymin <- m-sd(x)
    ymax <- m+sd(x)
    return(c(y=m,ymin=ymin,ymax=ymax))
  }

  p <- df %>%
    filter(!is.na(!! x)) %>%
    group_by(!! x) %>%
    ggplot(aes(!! x, !! y, fill = !! x)) +
    geom_violin() +
    scale_fill_manual(name = "Guest star", values = colors[3:4],
                      labels = c("FALSE" = "Playing an original character", "TRUE" = "Playing themselves")) +
    scale_x_discrete(labels = c("FALSE" = "", "TRUE" = "")) +
    labs(x = "", y = ylab, title = title)  +
    stat_summary(fun.data = data_summary, geom = "pointrange", color = "black",
                 show.legend = FALSE)

  if(!is.null(breaks)){
    p <- p +
      scale_y_continuous(labels = labels, breaks = breaks, limits = limits)
  }else{
    p <- p +
      scale_y_continuous(labels = labels, limits = limits)
  }
  p
}

plot.guest.ratings <- episodes %>%
  plot_violin(x = self, y = rating, ylab = "Rating score", title = "IMDb ratings",
              limits = c(1, 10), breaks = seq(0, 10, 2))

plot.guest.votes <- episodes %>%
  plot_violin(x = self, y = votes, ylab = "Number of votes", title = "IMDb votes",
              labels = comma, limits = c(0, 4000))

plot.guest.views <- episodes %>%
  plot_violin(x = self, y = us_views, ylab = "Number of US viewers", title = "TV views in the US",
              labels = unit_format(unit = "", scale = 1e+6, big.mark = ","),
              limits = c(0, 35))

ggarrange(plot.guest.ratings, plot.guest.votes, plot.guest.views, nrow = 1, common.legend = TRUE,
             legend="bottom", widths = c(0.95, 0.97, 1.08))